Solitons as the early stage of quasicondensate formation during evaporative cooling 
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We calculate the evaporative cooling dynamics of trapped one-dimensional Bose- Einstein conden- 
sates for parameters leading to a range of condensates and quasicondensates in the final equiUbrium 
state. We confirm that solitons are created during the evaporation process, but always eventually 
dissipate during thermalisation. The distance between solitons at the end of the evaporation ramp 
matches the coherence length in the final thermal state. Calculations were made using the classical 
fields method. They bridge the gap between the phase defect picture of the Kibble-Zurek mechanism 
and the long- wavelength phase fluctuations in the thermal state. 
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The quasi-one-dimensional (ID) Bose gas in elon- 
gated clouds of neutral ultra-cold atoms [l[ Q differs 
markedly from the Bose-Einstein condensate (BEC) in 
three-dimensional geometries. One of the most remark- 
able features is the presence of two characteristic tem- 
peratures when the trapped gas is cooled jH. Below Tc, 
the lowest mode becomes appreciably occupied Q but the 
phase coherence length cx 1/T is smaller than the size 
of the system. This is called a quasicondensate. Be- 
low a second temperature T^, grows to the size of the 
cloud, and the state is a true BEC. The phase coher- 
ence in these states have been extensively studied both 
experimentally d, 01 and theoretically^, |5H§|. 

In thermal equilibrium, the variance in phase and one- 
body density matrix have been calculated jH, [HI, [ol, [iq| . 
Their short-range behaviour gives the phase correlation 
length 1(f) near the center of the trap: 

piiz, z') = {¥iz)^iz')) « pi(0, 0)e-l---'l/'^ (1) 

In equilibrium, the phase in a single experimental reali- 
sation varies smoothly over length scales /^[iq|. 

On the other hand, phase fluctuations have also been 
predicted from the Kibble-Zurek mechanism flTl (KZM) 
after the onset of condensation, when the system is far 
from equilibrium. These fluctuations are seemingly dif- 
ferent in nature than those discussed above. During 
evaporative cooling, phase defects in the form of grey 
solitons appear when crossing the characteristic temper- 
ature Tc. They are born when local condensation occurs 
faster than distant regions can communicate to agree on 
a common phase. When the expanding initial phase do- 
mains meet, soliton defects form on the interfaces be- 
tween them. Therefore, during the formation of a con- 
densate the phase experiences sudden jumps at the tem- 
poral position of every soliton, and phase domains appear 
between them of a size equal to the separation between 
neighbouring solitons. A natural question arises whether 
these pre-formed domains are somehow related to the 



phase fluctuations in equilibrium. The aim of this paper 
is to show how the phase fluctuations in these two cases 
are connected. 

The number of solitons while crossing Tc has been pre- 
dicted as a function of the quench rate [6]. A calculation 
where chemical potential was quenched at T = demon- 
strated the KZM for a uniform gas [7]. Here we show that 
the Kibble-Zurek scaling also applies for a realistic model 
of evaporative cooling in a trap. We have simulated the 
non-equilibrium evaporative cooling dynamics to the sta- 
tionary thermal state using the classical fields method 
(CFM)[12]. The observations imply that the solitons are 
indeed the early stage of development of the final phase 
fluctuations at equilibrium. Let us now proceed to the 
details. 

We consider a single-species Bose gas in a trap. We 
use harmonic oscillator units with frequency ujz. The 
initial state of TV = 10^ atoms is generated in a canon- 
ical ensemble at a temperature of ksT = 360huj-.\l^. 
well above Tc in the final state. The description of the 
system in the CFM is in terms of an ensemble of clas- 
sical field amplitudes tlj{z) evolving under a generalised 
Gross-Pitaevskii mean field evolution equation from ther- 
mally randomised initial conditions. This is equivalent 
to a truncated Wigner description[l3|, but with quan- 
tum fluctuations omitted. It is a good approximation 
for the ID gas in our regime, where even the quasicon- 
densate phase fluctuations are dominated by the thermal 
component [2]. The evolution in time is given by 

idt^{z, t) = [H{z, t) - iV{z, t)] ^{z, t) , (2) 

with the Hamiltonian 

H{z, t) = -\j^ + ^(^' *) + 9iD\i^{z, |2 (3) 

including kinetic energy, contact interactions with 
strength gijj = 0.31 (corresponding to ^^Rb atoms in 
a 10 X 1000 X 1000 Hz trap), and the time-dependent 
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FIG. 1: The evaporative cooling potential V{z,t). Solid line: 
initially (t = 0), dashed line: at the end of the ramp (t = tr), 
dotted line: during thermal equilibration (t ^ tr). 
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external potential (shown in Fig. [T]) 
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V{z,t) = U(t) 
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The first stage {t < tr) is the evaporative cooling ramp. 
V{z) is a Gaussian dip of constant trap frequency {ujz = 1 
in the units chosen) near z = 0, and standard deviation 
^yU(t). The depth of the dip decreases linearly from 
Uo = 100 to Ur = Uo/3. The ramp time tr is varied to 
obtain different final states. 

After the ramp, {tr < t < tmax ) the evaporation is 
stopped and the gas is allowed to thermalise for a longer 
period to tmax = + 1000 in a much deeper potential 
(rising to a depth /7max = 10/7o) that becomes effectively 
harmonic in the region occupied by the gas cloud. It 
is also necessary to include losses in the wings of the 
potential: T{z,t) = F^o [V(z, t)/[/(t)]^ with F^o = 10, 
7 = 50. This loss acts as a high-energy knife in the 
region beyond about two (t-dependent) standard devia- 
tions of the Gaussian dip. This is necessary to realisti- 
cally model the experimental properties of the trap by 
preventing once evaporated atoms from returning back. 
A lattice of 1024 points on a length L = 120 is used. 
The number of trapped atoms we obtain for t > tr is 
always around N ^ 1300. The widths and shapes of the 
clouds agree very well with the radius ztf = 
chemical potential /i = [9(^id7V)^/32]'^^/^) predicted by 
the Thomas-Fermi approximation. 

Figure [2] shows the time-evolution of a single realisa- 
tion of the experiment in the quasicondensate {tr = 75) 
and BEG {tr = 400) regimes. Animations for tr = 250 
are in[l5|. We see the emergence of a great number of 
defects on the healing-length scale in the early part of 
the cooling phase (e.g. around t = 6 — 10 in Fig[2](c)). 
The density dips in Figs. [2]^a,c,d) display all the char- 
acteristics of grey solitons: passing through each other, 
turn-around near the edge of the cloud when the central 
density reaches zero, motion much slower than the speed 
of sound, phase jumps, and width in agreement with the 
healing length ^ = 1/ y/giD{po{x) - Pmin). Here po{x) is 
the local background density, pmin the minimum density. 
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FIG. 2: Density (a,c,d), and phase modulo 27r (b,e) for a 
single realisation during the cooling ramp. In (a-c) tr = 75 
leading to a quasicondensate as t ^ oo, in (d-e) tr = 400 
leading to a BEG. Top: at the end of evaporative cooling. 
Genter: The entire ramp time. Bottom: the loss of the last 
soliton, leading to BEG formation. 



The density profile of the dips is a good match to the 
solitonic solution (pmin + {po - Pmin) tanh^ [{z - ^o)/^]^ 
where the dip is centered at zq. Note the broadening of 
the solitons when they approach the edges of the cloud in 
Fig. [2](a,d). Very clear phase domains are seen between 
the solitons in Figs. [2](b), and account for the majority 
of the spatial phase variation. 

As cooling progresses, many of the solitons are lost at 
the edges, and the size of phase domains grows. We ob- 
serve that if the final state achieved is a quasicondensate 
(see Fig. [3]), some solitons remain at t^, while in the BEG 
regime (large tr) all are already gone by this stage. Fig. [3] 
quantifies the crossover from quasicondensate to BEG as 
a function of ramp time tr. 

The number of condensed atoms No{t) is calculated 
from the maximal eigenvalue of the one-particle density 
matrix pi{z^ z'). Here, this is calculated as (?/^*(z)?/^(z'))s, 
with (•)s denoting averaging over realisations within 
the initial canonical ensemble. The condensate frac- 
tion no = Nq/N, is shown as the 'x' data in Fig. [3l 
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FIG. 3: Fraction of atoms in the four lowest energy modes 
^0,1,2,3 (black), and phase domain width /1/2 (red) at the end 
of the evaporative cooling {t = tr) as a function of the ramp 
time tr. The horizontal line (red) shows the Thomas- Fermi 
width {2ztf — 2a/2/I) achieved as t y t-max 

). The ensemble 

consisted of 100 independent realisations. /1/2 > 2ztf occurs 
when phase coherence extends into the tails of the cloud. 

Phase correlations are characterised using the correla- 
tion function g'^^^z^z') = pi{z^ z')/^/pl{z^ z)pi{z\ z'). 
We use the width of the correlation function symmet- 
ric around zero to measure the phase domain size, using: 
g^^\—\li/2-) \h/2) — \ ^ Hence, comparing to 

©,^0 = ^1/2/ log 2. 

Fig. [3] shows the transition between quasicondensate 
and BEG behaviour [2]. In the quasicondensate, the phase 
domain width is much smaller than the cloud, and a 
sizable part, or even a majority of the atoms are not 
condensed into the lowest energy mode. In the con- 
densate, a single phase domain covers the entire cloud, 
and the condensate accounts for the vast majority of 
the atoms. Notably, the transition in phase correlation 
length /1/2 = ^</)log2 around tr = 200 is much sharper 
than in condensate fraction no- 

As the state is non-thermal during the evolution, as- 
signing a time-dependent temperature is moot, and we 
use the condensate fraction no (always a well defined 
quantity, unlike T) as the analogous parameter. For in- 
terpretation, several approximate results will be usefuljSl]: 
no ^ 1 — T/Tc, when no is not too close to zero, 
and where the characteristic temperature Tc is given by 
N = Tc log(2Tc). For TV = 1300, this is = 214. Hence, 
the temperature appearing in estimates of the coherence 
length jil, [HI can be replaced by 

T^T,(l-no). (5) 

The location of the BEG/quasicondensate transition in 
Fig. [3] can be compared using (j5j) with the prediction 
T0 = N/p^. The latter leads to nQ{T^) = l-N/{fiTc) ^ 
0.84, which compares favourably with the data. 

To verify that the Kibble-Zurek mechanism (KZM) 
is at work here, we wish to compare the scaling of the 
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FIG. 4: The mean number of solitons at tr (circles; error 
bars show standard deviation of the ensemble values). The 
solid line shows the scaling N oc 1/tr approximately expected 
for mean- field systems from Q , fitted to the data in the qua- 
sicondensate regime (tr < 200). The prefactor is ?^540. The 
inset shows a logarithmic scale. 

number of solitons (predicted in [6]) with the rate that 
the characteristic temperature Tc is crossed. We count 
solitons by fitting the density in the solitonic solution 
in the same manner as in 0, using a local background 
density. To make the comparison, we should also cal- 
culate the quench time of the relative temperature tq 
in which the KZM is expressed[6]. The prediction is 
iVsoiiton OC ^-(i+2-)/(i+-.) ^here f = -TJtq, and 
and uZr are critical exponents for healing length and re- 
laxation time at T = Tc, respectively. Here, because the 
system is certainly not in thermal equilibrium during the 
evaporative cooling phase, neither a value for Tc nor ^ 
is forthcoming. The most reasonable straightforward es- 
timation is that the quench time is proportional to the 
ramp time: tq oc tr- This leads us to use tr as the 
"quench time" parameter in Fig. |4] where the scaling is 
examined. 

We see that the KZM power law scaling holds quite 
well in the quasicondensate regime for power law expo- 
nent unity in 1/tr, which is the value expected for mean 
field [6]. As one might expect, the scaling law breaks down 
once the true condensate is reached for tr ^ 200, since 
the correlation length reaches or exceeds the size of the 
system in the BEG regime. This suppresses the slowing 
down near Tc that is an assumption of the KZM. 

In contrast to the elongated 3D case of no conden- 
sation front is apparent here. This may be an effect of 
the slow internal thermalisation of ID gases. We con- 
jecture that a front might still occur in ID if the local 
temperature were set externally by contact with an exter- 
nal reservoir, as e.g. in buffer gas cooling, but is unlikely 
in the type of scenario simulated here. 

Let us now consider the thermal equilibrium state at 
t = ^max- We observe that there are indeed no solitons 
present. The phase domain size calculated from g'^^^ is 



4 



quasicondensate: tr = 75 BEC: tr — 400 
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FIG. 5: Coherence length measures, (a-d): g^'^\z,—z) cor- 
relation function, (e): Phase domain widths as a function of 
condensate fraction no. Red dots: I1/2 values at equilibrium 
t = tmax; open squares: mean distance between solitons 
ai t = tr] dashed line: stationary state prediction for /1/2 
{h/2 = /0/log2 = {3/2ZTF log2)(iV/Te)(y^) from [2,]). 



gas is evaporatively cooled. 

In 2D, where vortices take the roles of solitons in ID, 
a similar sequence of defect creation and slow dissipa- 
tion into long-wavelength phase fluctuations should be 
active when crossing the BKT transition (17|. However, 
there, c-field simulations showed vortices still present in 
the thermal state(l8l.[l9|. while in ID the solitons always 
eventually dissipate. 

In conclusion, we have calculated the evaporative cool- 
ing dynamics of a ID gas through to the quasiconden- 
sate and BEC regimes in an experimentally realistic 
model. The simulations confirm the action of the Kibble- 
Zurek mechanism for soliton formation, and their sub- 
sequent decay into long- wavelength phase fluctuations. 
The match between inter-soliton distance at the end of 
the ramp and phase coherence length at long times shows 
that the solitons are the early stage of the formation of 
the thermal phase fluctuations. 
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this topic. This work was supported by Polish Gov- 
ernment Research Funds: N N202 104136, 2009-2011 
(EW,MG), N N202 128539, 2010-2012 (PD), N N202 
174239, 2010-2011 (KR). 



shown by solid circles Fig. [5je), where we also compare 
to the equilibrium theoretical prediction of [2] (dashed). 
The agreement is quite good deep in the quasicondensate 
regime, with a 2 — 3x discrepancy for 0.5 ^ no ^ 0.8. 
We note that experiments also tend to show coherence 
lengths longer than the theoretical estimate shown[3, lOj. 

Having calibrated our results with the two phase fluc- 
tuating situations, a connection between them can be 
made. We observe that the soliton creation timescale 
here is very short in comparison with the thermalisation 
time. A hypothesis is that initial contact of coherent 
phase domains formed via the KZM results in a soli- 
ton, which then converts into smooth phase fluctuations 
only on the long thermalisation timescale. Such a tem- 
poral mismatch between phase and density fluctuation 
behaviour is also reminiscent of recent experiments in ul- 
tracold Helium flG*]. A symptom of such a process would 
be to see the disappearance of solitons without change 
in the coherence length. Indeed, in the quasicondensate. 
Figs. [5](a,c) show no appreciable change from tr to tmax, 
despite the disappearance of solitons in the meantime. A 
"soliton" measure of the phase domain size is the mean 
distance d between them. If phase fluctuations are "born" 
as solitons, then d at earlier times should match the final 
phase domain width /1/2 in equilibrium. Fig.[5l compares 
d{tr) at the end of the ramp with at equilibrium. The 
match is remarkably good, strongly backing up the above 
hypothesis that solitons can be considered the precursors 
or seeds of the equilibrium phase fluctuations when the 
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